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We present a 3D quantum mechanical theory of radio-frequency outcoupled atom lasers from 
trapped atomic gases in the presence of the gravitational force. Predictions for the total outcoupling 
rate as a function of the radio-frequency and for the beam wave function are given. We establish 
a sum rule for the energy integrated outcoupling, which leads to a separate determination of the 
coupling strength between the atoms and the radiation field. 

For a non-interacting Bose-Einstein condensate analytic solutions are derived which are subse- 
quently extended to include the effects of atomic interactions. The interactions enhance interference 
effects in the beam profile and modify the outcoupling rate of the atom laser. We provide a complete 
quantum mechanical solution which is in line with experimental findings and allows to determine 
the validity of commonly used approximative methods. 

We also extend the formalism to a fermionic atom laser and analyze the effect of superfiuidity on 
the outcoupling of atoms. 

PACS numbers: 03.75.Pp, 03.75.Ss, 03.65.Nk, 05.30.Fk 
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I. INTRODUCTION 

The possibility of creating an atom laser analogous to 
an optical laser was considered immediately after the 
creation of atomic Bose-Einstein condensates (BECs). 
Atom lasers can be operated by continuously extract- 
ing small amounts of trapped atoms in a coherent way. 
Atom lasers provide an important tool to analyze the 
properties of trapped atoms and many experimental ap- 
plications are expected due to their coherence properties 
[l| . In principle they offer the possibility to monitor the 
evolution of a BEC without the need to switch off the 
trapping potential. 

The first experimental realization of a BEC output 
coupler was reported in Ref. |2|, where short radio- 
frequency (RF) pulses changed the hyperfine state of the 
atoms. The inhomogeneous magnetic trapping field sep- 
arated the atoms into trapped and outcoupled compo- 
nents. Using a series of RF pulses, a sequence of coherent 
atom waves was formed. 

A series of downward-falling output pulses analogous 
to a pulsed laser was demonstrated in Ref. |3| using an 
optical lattice. A BEC was loaded into a vertical stand- 
ing wave created with laser beams pointing in opposite 
directions. By lowering the depth of the lattice, phase- 
coherent atoms from different wells tunneled out of the 
traps and accelerated in the gravitational field. Simi- 
larly, it is possible to release an extended wave packet 
from a single optically trapped BEC by slowly lowering 
the trapping potential [; 
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A well-collimated quasi-continuous atom laser was 
achieved using a stimulated Raman transition as outcou- 
pling mechanism [6|. The two frequency Raman process 
imparts a momentum kick to the extracted atoms allow- 
ing directional output coupling. A sequence of overlap- 
ping matter wave packets were extracted from a BEC 
using repeated Raman pulses. Also the first continuous 
high flux Raman atom laser has recently been reported 



Using an extremely stable novel magnetic trap, Bloch 
et al. [3 demonstrated a quasi-continuous output cou- 
pler for magnetically trapped atoms. A weak RF field 
induces spin flips between trapped and untrapped hyper- 
fine states. The untrapped atoms fall in the gravitational 
field producing a coUimated atomic beam whose duration 
is determined by the condensate size. The possibility of 
continuous feeding the atomic source was demonstrated 
in Ref. Q . By using two different RFs on the same con- 
densate, the coherent spatial nature of the atom laser 
beam was shown in Ref. [l3|. The temporal coherence 
of atom lasers was investigated in Ref. [ly, and more 
recently also the second order temporal correlation func- 
tion [13. 

Since atom lasers generate coherent matter waves in 
the gravitational field of the earth, an accurate descrip- 
tion of the quantum-mechanical propagation in a linear 
force field is an im por tant part of a theoretical model of 
an atom laser [l3l [lj| . One dimensional models [13 are 
not sufficient for the characterization of the beam wave 
function, which requires a fully three-dimensional theory. 
Previous calculations of atom lasers rely mainly on nu- 
merical integration without including interactions [l6j, 
or employ semiclassical approximations for the propaga- 
tion in the presence of gravity and a mean field potential 

dill El. 



Recently, radial structures perpendicular to thegravi- 
tational field have been observed in atom lasers [la, llal ■ 
Similar structures have been predicted for smaller con- 
densates [Tj, I23, |21| , where they are linked to two-path 
interference in the presence of a linear force field. 

In this paper, we formulate and apply a theory of an 
atom laser, which employs the full 3D quantum mechan- 
ical propagator. After some basic definitions in Sec. ^ 
we review and extend in Sec. lIIII the analytic solution for 
the beam profile and the total current in the case of a 
non-interacting BEC. The analytic solvable model of an 
atom laser supplied by an "ideal" BEC forms the basis for 
the inclusion of mean field potentials in Sec. IIVI There, 
we derive a quantum mechanical multiscattering theory 
for the interacting BEC. We compare the quantum solu- 
tion to experiments and approximative methods like the 
Born series and the reflection approximation in order to 
judge the underlying assumptions and shortcomings of 
these approximations. The inclusion of a mean-field po- 
tential changes the emission rate of an atom laser, but 
the rate still obeys an important sum rule. The beam 
profile develops a radial substructure, which can be ob- 
tained from our quantum mechanical model and is in line 
with experimental observations. 

The effect of spatially anisotropic traps is discussed in 
Sec. where we show how one-dimensional and three- 
dimensional models are related. 

In Sec. I VII we extend the formalism to the current from 
higher modes in a harmonic trap. We apply it to a 
quasi one-dimensional Fermi gas and discuss the effect 
of fermionic superfluidity in both the current and the 
outcoupled density profile. 

The quantum source formalism provides a consistent 
framework for all presented calculations. A quantum the- 
ory for the beam wave function is of special importance 
for the coherent quantum control and tailoring of matter 
waves [l|. 



II. THE EMISSION RATE AND THE ATOMIC 
BEAM WAVE FUNCTION 

The output coupling of magnetically trapped atoms 
can be understood in terms of a spin fli p o f the mag- 
netic hyperfine quantum number mp |8|, llq (for *^Rb 
atoms, the F — 1 hyperfine level is commonly used). Ini- 
tially, the atoms in the state \F ~ l,mp = —1) are in 
an eigenstate of the atomic trap Hamiltonian in the pres- 
ence of a static magnetic field Bz and the gravitational 
field F = mg [g sa 9.81ms~^) along the z-axis 
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-mFgF^J-BBz - Fz 



(1) 



Here, gp denotes the Lande factor and /i^ the Bohr mag- 
neton. The inhomogeneous magnetic field of the atom 
trap is expanded in second order as a harmonic oscilla- 
tor potential for the mp = — 1 state. The gravitational 



field merely shifts the origin of this oscillator along the z 
direction and could be absorbed in the quadratic term. 
The application of an additional oscillating magnetic field 
with frequency i^ and amplitude B' adds to eq. JQ) the 
time-dependent potential V{t) = — ^iB' cos{vt) which 
causes transitions of the spin to the magnetic quantum 
number mp = (for simplicity we will not consider the 
m,p = 1 state). However, the \F = 1, m^p = 0) state is no 
longer an eigenstate of the trapping Hamiltonian, which 
only supports the mp = — 1 state. Instead its evolution 
is governed by the Hamiltonian 



grav 2yj^ ' 



(2) 



which leads to the propagation away from the trap and 
the formation of an atom laser beam. In the following we 
will assume a relatively weak amplitude B' of the oscillat- 
ing field, which makes it possible to deplete the ground 
state of the trap over some time. Since our main interest 
is the determination of the atom laser wave function and 
the emission rate as a solution of the Schrodinger equa- 
tion, we will treat the radio wave as a classical radiation 
field. 

At this point we still have to solve for the time- 
dependent eigenstates ip[,.^p and V'grav of a coupled sys- 
tem, which have an energy difference of AE = i?grav ^ 
i^tiap and are coupled via 7 = \iB' 

{ihdt - iJgrav)V'grav(r, i) = 7 C-'^^*/VU(r' *) (3) 
(ihdt ~ iftrap)V';rap(r, t) = J e+'^^*/ Vgrav (r, t) . (4) 

Here, we employed the rotating wave approximation. We 
split off the time dependence of the states 

(5) 
(6) 

(7) 
(8) 



^grav(r,t) = e-'^--*/Vgrav(r) 

Cap(r,t) =e-'^-p*/VU(r)> 
in order to obtain the stationary equations 



i,^grav -'^ grav j'^ grav \^ ) 
{Etre,p - i?trap)V'trap(l") 



= 7 V'trap(r) 
7 V'grav(r). 



Now we will use the assumption of a weak coupling in or- 
der to replace the state ^t^ap (i") by an eigenstate of eq. (Q) 
denoted by V'o(r). Thus we break the coupling between 
the two equations and are left with the evaluation of the 
stationary Schrodinger equation 

(-Bgrav - -ffgrav)V'grav(r) = 7 V'o(r) (9) 

in the presence of an inhomogeneous source term 



cr(r) =7V;o(r). 



(10) 



The time-independent source term is akin to the steady 
state solution for the atom laser beam after the damp- 
ing of transient effects due to the initial switching on. 
In the following, we restrict the discussion to the sta- 
tionary case. The inhomogeneous equation is readily 



solved by using the energy-dependent Green function 
Ggrav(r,r';£') for -ffgrav The Green function is the solu- 
tion of the Schrodinger equation for a point inhomogene- 

ity 



{E - i?grav)Ggrav(r, r'; E) = Sir - r'), 



(11) 



and thus we obtain the wave function emitted from an 
extended source cr(r) by a convolution integral: 



^; 



grav 



{v-E) = j dv'G^,^^{v,v';E)a{v'). 



(12) 



In the following we will suppress the subscript grav, since 
we are always interested in the properties of the atomic 
beam. We take the energy E as the difference between 
the energy of the radiation field hv minus the Zecman 
splitting between the levels mp = — 1 and mp ~ 0: 



hv - {E^ 



Eu2ip ) ■ 



(13) 



Notice that due to continuous spectrum of i/™!^^ , an 
output coupling is possible for a continuous range of en- 
ergies. However, for i<^ — > ±oo we will see that a sum 
rule for the energy integrated outcoupling rate enforces 
a vanishing outcoupling rate. We now proceed to calcu- 
late the total current J{E), which denotes the number 
of atoms released per second. We do this by defining a 
current density associated with the wave function H12|) 

j{r;E) = ^^[^bir;ErV^iv;E)]. (14) 

By using eq. ^ it is straightforward to derive the equa- 
tion of continuity for this stationary problem: 



V-j(r) = --3{a(r)*V(r)}. 



(15) 



As wanted, the inhomogeneity models a constantly emit- 
ting particle source. The integration over a surface en- 
closing the source yields a bilinear expression for the total 
probability current: 

J{E) ^ ^'^^U Av j Av' a{vr G{v,v'-E)a{v'] 

(16) 
Usin g th e alternative representation of the Green func- 
tion ll El, 



G{v,v';E) = {v 



1 



E- H 

we can rewrite eq. ((TC)l as [Ij, |23 
27r 



- m5{E - H) 



J{E) 



fi 



(t\5{E - H)\(j) 



(17) 



(18) 



The quantity {y\6{E — H)\y) is the local density of states 
(LDOS) of the Hamiltonian H . For initial states normal- 
ized to (-00 1 ^o) ~ N an important sum rule [l^l follows 
from eq. ((TH|l : 



/■°° 27r 

j^^AEJ{E) = —{a\a) 



2ttj^N 



(19) 



One consequence of the sum rule is, that for any Hamil- 
tonian we can determine the interaction strength by sim- 
ply summing up the total current at different energies. 
In Ref. [13 this was used to check the experimental re- 
ported coupling strength 7. Also the finite value of the 
integral in eq. (|19(l restricts the output coupling to a spe- 
cific energy range. In the next sections, we derive analytic 
expressions for the Green function G{r,r';E). 

It is worth mentioning that the source formalism in- 
troduced here is analogous to the first order perturbation 
theory in the coupling Hamiltonian used for fermionic RF 
outcoupling in 0ISI2I. 



III. THE IDEAL ATOM LASER 

In this section wc discuss an extension of the theory 
for an ideal atom laser presented in Ref. jTj|. The ideal 
case forms the basis for the inclusion of interactions into 
the theory in Sect. lIVI For the Hamiltonian of the linear 
force field (O there exists a closed analytic expression of 
the Green function EHMISi: 



>-rgravVr,r ,_C/j — ^^2 I ,1 : 



2r 



(20) 



where 



u± = -p[2E + F{z + z')±F\r-r'\], (21) 

(22) 



/3 = [m/{Ah'^F'^)Y^\ 



and Ci(x) = Bi(x-) + iAi(a;). The application of this 
Green function to an ideal atom laser from the ground 
state of a BEC is discussed in detail in [l3l a nd extended 
to include vortices in rotating BECs in |Tj|- In Sec. 4.2 
of |l3 | , an analytic solution for an atom laser originating 
from a non-interacting isotropic BEC with N atoms of 
Gaussian form 

a(r) = 7ViVa-3/2;r-3/4c--V2a^ ^ (23) 

is derived. Using scaled variables 

{^,vX)=f3F{x,y,z), e = -2pE, a^^Fa, (24) 
C = C + 2a^ p^=C^+t^^ + C^ ^=e + 4a^ (25) 

and the special functions defined in [lj|, App. B 

Qi(p,C;e) = -:^[Ci(e-C-p)Ai'(e-C + p) 

^Ai(e-C + p)Ci'(e-C-p)] (26) 



Qr"(p,C;e) = - 



du- 



-1(^-0+^ 



27r3/2y3/2 



Qii(6) = [A:^{e)f-e[ki{e)f 



(27) 
(28) 



one obtains analytical expressions for the total current 
([a, eq. (73)) 

Jideai(i^) = I (i{pFfK{efq\,{i), (29) 
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Detuning frequency Au [kHz] 

FIG. 1: Total current (per atom) as a function of radio fre- 
quency detuning Av = E/h from eq. 1291 (for the sign of 
the detuning see the endnote [3011'). Also shown is the reflec- 
tion approximation (dashed line). The area underneath the 
curves is fixed by the sum rule H19|l . Parameter: a — 0.4 /xm, 
j/h = 100 Hz, F — mjn^g, with g — 9.81 m/s^, and 
niKh = 87 u. 



described by an Gaussian profile [13, right panel 
of Fig. El 

We are not aware of the experimental observations of 
atom laser beams from small condensates, although high- 
quality magnetic microtraps are able to produce the re- 
quired BECs [33. 

For larger condensates one can derive approximative 
expressions for the total current. Of special simplicity is 
the so-called reflection approximation, which was devel- 
oped in the theory of Franck-Condon factors ^2I- I* is 
analogous to the local density approximation commonly 
used in BEC theory and consists of neglecting the ki- 
netic energy term in the Hamiltonian 101. Now eq. H18() 
simplifies to 



'^ideali^l 



^'^ ^dr|cr(r)|^(5(£;- 



Fz) 



(32) 



and the beam wave function 

V'idcai(r; E) = -4/3(/3F)3A(e) [Qi(p, C; e) + Q^'Xp, C; e)] ■ 

1-1 ^^°^ 

In addition to eq. (72) in Ref. jlj|, we explicitly add the 

near field contribution denoted by Q"°'''. For the non- 
interacting case, this term has no influence on the beam 
profile outside the condensate region, whereas we have to 
include it for the interacting model in the next section. 
In both expressions, the source strength A(e) is strongly 
energy- and size-dependent: 



A(?) = VA^7(2V^a)'/'e2"'('^-4"V3) 



(31) 



Let us briefiy summarize the main features of the ideal 
atom laser model (see also Section liVD|) : 

• As shown in Il3 |. Sec. 3.3, the beam wave function 
can be mapped back to a virtual tunneling point 
source, which is closely related to the Green func- 
tion. 

• For small condensates (radius a < "]fj^ in the di- 
rection of the gravitational field, i.e. about 0.5 /im 
for a Rb BEG), additional modulations in the total 
current appear and the beam wave function devel- 
ops an interference structure. The emerging pat- 
tern can be explained in terms of two-path inter- 
ference in the linear force field. A typical sequence 
of beam profiles for different detuning energies is 
shown in left panel of Fig. [3 



For large condensates (radius a > ™y^° in the di- 
rection of the gravitational field), the energy de- 
pendency of the total current reflects the density 
distribution of the source. The beam wave function 
is featureless over the whole energy range and well 



The reflection approximation states that the current is 
proportional to a slice through the initial density distri- 
bution along a plane height z = E/F. It is possible to 
justify this approximation as a limit of the quantum solu- 
tion (see [13, eq. (40)). In principle one can also calculate 
quantum corrections to eq. (|S^ (see [23, I^J, 113, Hjl ) , but 
the resulting (asymptotic) series can diverge for a finite 
number of terms (see Ref. |33, Sec. VI). 

The range of validity is limited by the requirement for 
a large spatial overlap between the initial wave function 
and the outgoing wave function in order to average over 
the oscillations of the Airy function in Eq. H20|l . If the 
width of the initial wave function is smaller than the first 
oscillation period of the Airy function (given by approx. 
1/{(3F)), the oscillations in the outgoing wave function 
carry over to the total current (see Sect. lVJI . 

In Fig. n we compare the reflection approximation with 
the quantum mechanical result eq. H29() . 

While the reflection approximation can in certain lim- 
its reproduce the total current distribution, it cannot 
yield information about the atomic beam profile. The 
slicing picture may be suggestive for the idea that atoms 
are only outcoupled at the slice given by the condition 
z = E/F. However, this assumption is not supported 
by the quantum mechanical result (|30|l . Gontrary to the 
picture of atoms leaving the condensate with zero mo- 
mentum along a slice (which would in a (semi-)classical 
picture imply that the radial profile of the beam at all 
distances from the BEG is identical to the density profile 
of the BEG), the beam profile spreads out as shown in 
Fig. 121 The summation over the infinitely many starting 
points distributed over the complete BEG is represented 
by a single point above the center of the condensate, and 
not by a planar surface at z = E/F. The semiclassical 
picture is further discussed in Sec. IIVDI 
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FIG. 2; Atom laser beam profile at different detuning energies 
t^v = Ejh from eq. 1301 for a relatively small BEG (left panel, 
a — 0.4 /xm) and a larger BEG (right panel, a = 0.8 /xm). For 
the sign of the detuning see the endnote l30ll . Shown is a 
cut through the middle of the beam along the vertical axis 
from 1 mm to 2 mm below the BEG for the detuning fre- 
quencies (0, 1, 2) kHz. The beam profile widens and develops 
a transverse substructure for larger energies. There is rota- 
tional symmetry about the vertical-middle axis of each profile. 
Parameter: a = 0.4 /im, F — niRb g, with g — 9.81 m/s^, and 
mRb = 87 u. 



Simultaneous outcoupling with two different radio 
frequencies. 

Experiments which outcouplc atonilascrs from a BEC 
with two different radio frequencies at the same time |inl | 
have shown the appearance of longitudinal interference 
structures. Using the quantum mechanical theory of the 
previous section, the resulting atom laser beam is de- 
scribed by the coherent superposition of two stationary 
beams with different energy originating from the same 
virtual point source: 



IV'rrai(r,t; ^1,^2)1 = h/'idcai(r;/iz^i)e 



ihi'it 



(33) 



The superposition of the two beam wave functions leads 
to a time-dependent oscillation of the density profile, 
which reproduces the observed longitudinal interference 
structure as shown in Fig 01 The use of multiple radio 
frequencies provides an important tool for tailoring the 
atomic beam wave function. 



IV. MEAN FIELD EFFECTS IN THE ATOM 
LASER 

In this section we extend the quantum theory to in- 
clude interactions between the BEC and the emitted 
atom laser beam (interactions within the atomic beam 
can be neglected since the density is much smaller than 
inside the BEC). A commonly used approach to include 
interactions in the description of a BEC is the addi- 
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FIG. 3: Squared beam wave function for simultaneous output 
coupling with two different radio frequencies. Shown is the 
density given by eq. 1331 . The outcoupling frequencies Ai/1,2 
are (a) ±0.5 kHz, (b) ±1.0 kHz, and (c) ±2.0 kHz. The 
number of longitudinal interference fringes is proportional to 
the difference in the detuning frequencies. Parameter: a = 
0.8 ^m, F = mB.bg, with g — 9.81 m/s^, and rriBb ~ 87 u. 



tion of the (repulsive) mean-field potential via the Gross- 
Pitaevskii equation 



Vgp(i-) =gsc|i/'o(r)|' 



(34) 



where gsc denotes the interaction strength and is related 
to the scattering length Qsc and the number of atoms in 
the BEC N via 



5sc = Anh^Nasc/m. 



(35) 



In principle the density dependent term leads to a non- 
linear Schrodinger equation. However, for the theory of 
the atom laser we will treat the BEC (and |?/'o(r)| ) as 
unchanged during the outcoupling process. Therefore we 
just have to modify the Hamiltonian for the propagating 
state by the mean field potential of the BEC: 



Zm 



(36) 



The additional repulsion will lead to a broadening of 
the beam compared to the non-interacting case, as we 
will show next. The repulsive interaction also leads to 
a change in the BEC density distribution itself. In the 
following we will retain the Gaussian form of the conden- 
sate, but the half-width a of the BEC should be viewed 
as a parameter, which can be obtained i.e. by minimiz- 
ing the Gross-Pitaevskii energy- functional [33 . For large 
condensates, a Thomas-Fermi like profile of the density 
is more appropriate than the Gaussian approximation, 



whereas for condensates in strongly confining traps a 
Gaussian profile is a fairly good approximation [31|, be- 
cause larger trapping frequencies increase the ratio of the 
kinetic energy vs. the interaction energy for a constant 
maximum condensate density (see Ref. |33, Sect. 6.2). 
However, the following discussion and comparison of dif- 
ferent methods for the atom laser rate and profile is in 
principle not limited to an initially Gaussian density dis- 
tribution. 



A. Quantum Theory 

The Green function of the Hamiltonian H3()|) is not 
available in analytic form. In principle, the Born series 
could be used to construct the Green function: 



90 



Ggp 



G 



grav 



^Siav^P*-^gra 



Ggrav VgP Ggrav VqP G 



^grav 



(37) 



As we will see below, for a typical BEG this series con- 
verges very slowly. An alternative approach consists in 
decomposing the mean field potential in terms of a 5- 
lattice 



JV 



F,(r) = ^K3p(r,)Ar<5(r-r,), 



(38) 



where Ar denotes the volume element of each lattice site. 
In order to mimic a continuous potential, the lattice spac- 
ing must be smaller than the typical oscillation length of 
the Green function, which is given by A « \/{(3F). Nu- 
merically, convergence has been checked by a set of cal- 
culations with subsequently reduced lattice spacing. For 
our calculations, we used a spacing of 0.15/im, which is 
about A/4. The (5-lattice is algebraically solvable via the 
transition (T-) matrix method (for a compact derivation 
see Ref. [S^, App. D). The T-matrix involves only the 
known Green function of the linear potential 



T{E)-^ = 



^trgrav(rj, ffc; E) 

[FGp(r,)Ar]- ' 



^grav \^ 3 



{r„E) j^k 



(39) 



and the renormalized Green function Gg°™(r; E) for r 
r' (la, (D21)): 



GgravXl":-^) 



m/3F 



[uCi{u)Ai{u) - Gi'(M)Ai'(u)] 



with u 
reads 



(40) 
-2(3{E + Fz). The resulting Green function 



GGp(r,r';S) = Ggrav(r, r'; £;) 



(41) 



N 



^ Y, Gg,.av(r,r,;S)T,fc(S,r,,rfc)Gg,av(rfc,r';S). 
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FIG. 4: Total current for a BEC of iV = 500 atoms as a func- 
tion of the detuning frequency IS.V = E/{2Tifi). We compare 
the results using the T-matrix method 11421 (with 1309 lattice 
points spaced at distances of 0.15/^m), the reflection approx- 
imation (1461 . the first order Born approximation 1451 . and 
the non-interacting result 1291 . Parameter: N = 500 atoms, 
tSscat = 5.77 nm, a = 0.8 /^m, 7/ft = 100 Hz, F — mab g, with 
g = 9.81 m/s^, and rriRb = 87 u. 



Using the Green function Ggp, wc proceed to calculate 
the total current from cq. (|16|l : 

Jgp(^) = -|3[(a|Ggravk) 

+ Y^ Tjfe(^)V'idcai(rj; £;)^idcai(rfe; -B)] 

= JidcAE) + J^p{E). (42) 

The sum rule p9|) enforces that the changes of the current 
vanish upon integration over the outcoupling frequency: 



/oo 
dEJip{E)=Q. 
-00 



(43) 



Similarly, the beam wave function (see eq. (|12|) 'l is given 
by the sum of the ideal profile and an interaction term 

ipGpir]E) = Vidoai(r; -E) 

+ Y Tjfe (-B)Ggrav (r , r, ; -B)^idcai (r^ ; E)] 



j.k 



V'idcai(r;£;)+V'Sp(i';^)- 



(44) 



B. First order Born and reflection approximation 

Within the T-matrix approach, we obtain the first or- 
der Born approximation by setting 



Ljfc 



™(ii;)=,5,feK3p(r,)Ar. 



(45) 



In general, the first order approximation is not sufficient 
for an accurate description, as shown in Fig. ^ 



In analogy to the non-interacting case, we can include 
the mean field potential in the reflection approximation. 
The previously planar slices arc now distorted, depending 
on the density of the condensate wave function: 



kHz 1 kHz 2 kHz kHz 1 kHz 2 kHz 

-1mm— I ^ I I -1mm — 



JhfiE) 



277 r 
'nj 



dr\a{r)\^d{E + Fz-VGp{r)). (46) 



C. Comparison of the quantum theory and 
approximative methods 

In Fig. 2] we compare the results of the different meth- 
ods for the total current as a function of the detuning 
frequency. We choose a Gaussian condensate with half 
width a = 0.8 fira, for which we expect a tunneling be- 
havior up to detuning energies of Fzq ~ 8 kHz. The 
ideal (non-interacting) total current reflects the Gaussian 
density profile of the condensate and attains therefore a 
symmetric shape. The inclusion of the mean-field poten- 
tial shifts the maximum of the total current to higher 
energies, as shown by the quantum mechanical T-matrix 
calculation. The reflection approximation works surpris- 
ingly well, whereas the first order Born approximation 
gives a misleading result with an additional hump. No- 
tice that the area underneath all curves is the same, as 
required by the sum rule H19|) . The shift to higher values 
of the detuning frequency in the maximum of the out- 
put couplingis in agreement with experimental results 
reported in [3, Il3 (note that the definition of the sign 
of the detuning frequency used here is opposite from the 
one used in the experimental work). 

Due to the effective negative initial kinetic energy, 
the beam profile of the non-interacting atom laser has 
a Gaussian profile, without a radial substructure (see 
[13). However, the presence of the repulsive mean field 
potential affects the beam profile considerably as shown 
in Fig. [3 In the transverse direction a substructure de- 
velops, which has been observed experimentally p!a.ll9l|. 
In a simple one-dimensional picture, the widening of the 
beam profile as compared to the non-interacting case has 
been attributed to the repulsive hump in the potential 
acting as a diverging lens |l7l |. The three-dimensional 
quantum mechanical picture is more involved, since the 
T-matrix approach includes multiple scattering events 
and no simple semiclassical interpretation in terms of tra- 
jectories is available. 



D. Semiclassical models for the beam profile 

The appearance of an interference structure in matter 
wave experiments can be linked to the possibility of mul- 
tiple paths from the source (or emitter) of the wave to the 
location of the detector. In the presence of a linear force 
field, the classical double slit experiment for electrons |33 
was carried out by Blondel et al [i^, |41| without actually 
constructing a material double slit. The uniform field 



-2 mm- 



Hi 



I I 

0.1 mm 



-2 mm- 



I I 

0.1 mm 



FIG. 5: Left panel: Atom laser beam profile from Eq. 14411 at 
different detuning energies Av = E/h including interactions 
(A'^ = 500 atoms, ascat = 5.77 nm), other parameters as in 
Fig. 2] The beam profile broadens and develops a transverse 
substructure for larger energies. Right panel: Atom laser pro- 
files for the same detuning frequencies and condensate width, 
but without interactions. The transverse interference pattern 
is not present. For the sign of the detuning see the endnote 
1301 . Shown is a cut through the middle of the beam along the 
vertical axis from 1 mm to 2 mm below the BEG for the de- 
tuning frequencies (0, 1, 2) kHz. There is rotational symmetry 
about the vertical-middle axis of each profile. 



environment provides (for positive initial kinetic energy) 
a region, in which there are two paths connecting the 
source with the target [43 . In the semiclassical approxi- 
mation of the energy Green function all classical allowed 
paths carry a complex amplitude (determined by the clas- 
sical action) and are added coherently [l^ 143. 144. l45l | . 

In contrast to the classical analysis of the trajectories 
from a single point in space to another point, a spatially 
extended source region, like a BEG, seems to require the 
addition of infinitely many paths leading from every point 
of the source region to the target point. Remarkably, 
the single point interference pattern is not destroyed by 
this averaging process. The reason is that, similar to 
the technique of virtual point sources in optics, one may 
replace the extended BEG by a single point source which 
is located at a distance 



raFa^ 



zq 



9 



(47) 



2r 

above the center of the BEG [13 (we used a ~ 
^J h / {mio)) . The focal point of the parabola given 
by the trapping potential (converted to spatial units) 
VtTap{z)/F = jpmuj'^z'^ coincides with the location of 
the virtual point source at a distance zq from the con- 
densate. Due to the shift upwards in the gravitational 
field, the initially available kinetic energy is reduced by 
the potential energy at the shifted location 



Skin = hAiy - 1/gp(0, 0, zq) - |i^zo| 



(48) 



Here, we also included the potential term due to the mean 
field of the condensate atoms, which creates a hump in 



the otherwise planar potential surface of the linear grav- 
itational field. For a large condensate, the initial shift 
I zq I 3> a leads to a virtual point source which is actually 
located several half widths a apart from the center of the 
BEC. The resulting initial kinetic energy is negative for 
detuning energies in the range of energies for which we 
expect a significant total current /lAi/ < Fa. For the 
semiclassical analysis we note that no classical trajecto- 
ries exist up to the turning surface, which is given by the 
implicit equation 



/lAi^- Vgp(0,0,2)- \Fz\ = 0. 



(49) 



In principle, one can study the classical trajectories which 
start from this caustic surface and end at a given target 
point. The possibility of multiple trajectories leads to 
a coherent sum over the corresponding classical actions 
[l7lll9| . However, a caustic surface is not an ideal start- 
ing point for a semiclassical analysis, since the specifica- 
tion of an initial position and simultaneously a definite 
momentum p = is not compatible with the uncertainty 
relation J43. The unknown initial phase and weight of 
the manifold of classical trajectories starting from the 
caustic surface presents another difficulty for a well de- 
fined semiclassical description. 



V. GEOMETRY EFFECTS IN THE ATOM 
LASER 

So far we have only considered spherically symmet- 
ric clouds of atoms. In principle one can tunc the mag- 
netic trapping potential and thus vary the frequencies of 
the harmonic trap. Experimentally it is possible to ob- 
tain quasi one-dimensional (ID) systems |46ll47l | by tun- 
ing one of the frequencies of the three-dimensional (3D) 
trap to a value which is much smaller than the other fre- 
quencies or by using optical lattices to create arrays of 
smaller ID systems [i^. A ID gas with the long axis 
in the direction of the gravitational field is characterized 
by cjx <C w^ = {uj^. + cjy)^/^. The propagation occurs in 
the three dimensional space and thus it is not sufficient to 
consider merely the propagation in the gravitational field 
along the beam axis. While the total energy is conserved, 
it falls into parts related to propagation in the direction 
of the gravitational field and in the perpendicular plane 
respectively. In the following we analyze the relationship 
between three-dimensional and one-dimensional calcula- 
tions. For simplicity we do not include mean field effects 
in the calculations. 

We consider the current and outcoupled wave function 
for an initial state of the form 



?/'i„i(x,2/,z) = (r|i/'ini> = i'Q{x,y)il)f){z), 



(50) 



where tpo^z) and 'ipQ{x,y) = Xj sfna^e^^^ +^ )/2ai are 
the ground s tates of the harmonic oscillator with half 
width a — ^jTi/ [mujz) and a±^ ~ y^h/{'muj±). The con- 
nection between the ID and the 3D Green function of the 
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FIG. 6: Total current (per atom) as a function of the detuning 
frequency Av — E/h. Shown is the 3D total current H53|l for 
several geometries. Parameters: ^/h = 100 Hz, rriNa = 23 u. 



linear gravitational field is given by a convolution integral 
over the transverse momenta k and k' 

Gl?^,{v,v'-E) = ^ j j dkdk' e-'^^^-^'^-'^'^y-y'^ y. 
Gl?^Jz,z';E- ^'^^l^_^'"A , (51) 



2m 



where (|3|, eq. (B2)) 



GlEjz,z';E) = -47r/3^FCi(u+)Ai(u_). 



(52) 



For the total current (|16|l we evaluate the expectation 
value of the Green function with respect to the initial 
state ■i/'ini- The resulting convolution integral reads 



J-^^(S) = 2ai 



d/cifc_Le""i''ij^^ ( E- 



n^kl 

2m 



(53) 



where 



1^,^, _ 16(3TT^/'^a-f^N _i6„6 



J^^ (E) 



a°+4a^e , 



AiCeY (54) 



is derived in eq. (|65|) . 

In Fig. El we show the 3D total currents for a big 
condensate (|zo| > a or Fa > 2hL0z, see eq. (|T7J) ') 
whose width along the gravitational axis is a ~ 2.1 /im 
{ujz = 27rl00IIz) and for a small BEC with a = 0.66 /im 
{u!z ~ 27rlkIIz). Fig. shows the current for Na atoms 
instead of Rb, since then the condition for small conden- 
sates is fulfilled for snraller trapping frequencies. In the 
reflection approximation eq. (|32|) we expect that the out- 
coupling window is just determined by the condensate 
extension along the gravitational field and not changed 
for the one-dimensional and three-dimensional case. 

However, the convolution integral (|53|) of the ID 
current with an exponentially decaying function whose 



width is proportional to uj± predicts in general a differ- 
ent form of the ID and 3D current. This effect can be 
clearly seen in Fig. El Only for a large condensate the 3D 
current retains a Gaussian shape (see the dotted line in 
Fig. ^ . A tight confinement in the transverse direction 
results in a spread of the total current as shown by the 
dot-dashed line in Fig. |H1 

For small condensates the situation changes com- 
pletely. Now the total current in ID is modulated by the 
zeroes of the Airy function in eq. (|54|l . The convolution 
with the transverse momentum tends to wash out the ze- 
roes if a;_L is big enough (dashed line in Fig.|SJ), but there 
might still exist regions of almost complete suppression of 
the total current within the outcoupling energy window 
(see the solid line in Fig.EJ. 

One can also separate the expressions for the ID and 
3D outcoupled wavefunctions. The outgoing wave func- 
tion perpendicular to the gravitational force is given by 



hyperfinc states. Fermionic superfluidity relies on the for- 
mation of pairs of attractive atoms. In the limit of weak 
interactions the system can be described by the Bardeen- 
Cooper-Schrieffer (BCS) theory developed in supercon- 
ductivity that relates the order parameter to the bind- 
ing energy of the paired atoms (gap). A beam of atoms 
coherently outcoupled from a trapped gas preserves the 
properties of the initial state of the atoms. We explore 
the effect of quantum degeneracy and fermionic super- 
fluidity in the outcoupled beam density profile and the 
total current. For the sake of simplicity we consider quasi 
one-dimensional systems. Furthermore, it is reasonable 
to expect that the effect of superfluidity in the outcou- 
pling increases when the superfluid gap is highest along 
the direction of gravity. 



A. Excited modes falling in the gravitational field 



i^onti^^y) 



1 



dfcdfc' / / dx'dy' X 



4a 1 7r5/2 



/2 , /2 

-ik{x-x )-ifc {y-y ) V 



(55) 



= / dk^^-J-c 2^ 3o{k±\/x'^ +y^), 



where Jq denotes the Bessel function of 0th order. The 
corresponding 3D wave-function reads: 



V'oTt(r;i?) = 



a±_ 



l-K 



dfc^ fcj_e 2 Jo(fcj_-\/aJ2+y^) X 







f-.(=;i=-^) 



(56) 



In this section we present a one dimensional calcula- 
tion of an excited state of a harmonic trap falling in the 
gravitational field. The current is calculated from the 
Franck-Condon factors of the initial and final eigenfunc- 
tions. The initial Hamiltonian of a ID harmonic oscilla- 
tor reads 



f? d2 1 



(58) 



whose eigenstates and eigenenergies arc given by 

Er^ = a^j^n-F-), (60) 



where 



Ai(?)Ci(e-2C~). 



Q°+2Q^e 



(57) 



Similar to the total current, the three-dimensional 
wavefunction is given by the convolution of the one- 
dimensional wavefunction and a exponentially decaying 
function whose width is proportional to uj^. 



VI. CURRENT FROM HIGHER MODES IN A 
HARMONIC TRAP. FERMIONIC ATOM LASER. 

In this section we consider the current and the out- 
coupled wave function from excited modes in a harmonic 
trap. In particular we consider a fermionic gas at T = 0, 
but the formalism presented here can also be used to 
analyze the contribution to the current and outcoupled 
beam of the higher modes in a trapped BEG. 

Quantum degeneracy was demonstrated for a trapped 
cloud of fermionic alkali atoms in [i^l • A gas of ultracold 
fermionic atoms becomes superfluid [SOl Isil 153 . 153. |5J| 
by tuning the interaction strength between two different 



where a = ^Jhjrriujz- The final state Hamiltonian is given 
by the ID version of iJgrav defined in Eq. (0) 



--4^- 



(61) 



The Hamiltonian Hi has a continuous spectrum and the 
eigenfunctions can be labelled by the energy |55| 



^'e{z) = (2|V^) = 2VFpAi U2(3F{z + j)\ 



(62) 



Inserting the complete set of continuum eigenfunctions 
in eq. (|18|l yields the following expression for the total 
current originating from an initial state j^'n) 



4"(i?) = 



27r72 



d^g,av<5(-B - ^g,av)|(^„|^£)r. (63) 



The needed overlap integrals (which are Franck Gondon 
factors) are conveniently calculated by adapting a recur- 
sive method developed in Ref. [^ (see also Add. IX)| : 



{i^n\^E) 



8v^/3a /lea*^ .^i;. 2 , 
cxp — AI3Ea X 



2"n! ""'' \ 3 
A'(r^,2,-8a^-2/3£;-H4a^-4a). (64) 
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FIG. 7: Total 3D current for a quasi-lD Fermi gas (eq. (|75J) 
as a function of the detuning frequency ^v — E/h. Pa- 
rameters: N = 70 spin-polarized atoms, ujz = 27rl.2 kHz, 
uj± = 27rl20 kHz, Ep/h = 84 kHz, T = 0, 7//1 = 100 Hz, 
rn,K ~ 40 u. 



Here, iir(n, a^, a^,7L, (5l) is given by eq. (|A5|) . The ID 
current becomes 



j,^.^(i?) = ^i(vj^i.)p. 



(65) 



The outcouplcd wave function can be calculated from 
eq. (|nj) 



^l^,{n,z;E) = jJ dz'GMz,z';E)Mz')- (66) 

As shown in the appendix, a closed expression is given 
by eq. !(K^ 



4>ti:tin.z;E) - B{n,z,E)x 



(67) 



K{n, 2, -8a^ -2/3£; + 4a^ -4a), 
where B{n, z, E) is defined in Add. 1X1 

B. Long axis in the direction of gravity 

We consider the current and outcouplcd function for 
an initial state of the form 



V'ini(a;,y, 2) = (rl^/w) = i^a{x)^Q{y)^n{z) 



(68) 



where V'i are the eigenstates of the harmonic oscillator. 
Thus we consider an arbitrary excited state along the 
z direction while the wavefunctions are gaussians in the 
perpendicular plane. The total current is given by a con- 
volution integral with the ID current of eq. H65() : 



/•OO 

Jl"{E) = 2ai / dfcj 
Jo 



k 



±e 



-J, 



ID 



E^ 



n^kj 

2m 



The corresponding 3D outcouplcd wave-function is ob- 
tained by using Eq. ifS^ 



V'oTt(r;i?) 





L 


dkj_ 


fc^e-°"=" 


Jo(fc±v^ 


+ 


y2)x 


i^l^. 


[n,z 


■E- 


2m '' 






(70) 




C. 


Fermi gas 









(69) 



The wave function of a Fermi gas at T = is a Slater 
determinant of the product state of the wave functions of 
N atoms. The density profile and the current are just a 
sum of the currents and wave functions of the individual 
states of the atoms 

nout(r; E) ^Y. /(^")lV'out(n, r;E + S„)|2 (71) 

n 

J(S) = ^/(S„)J„(S + S„), (72) 



where /{E^) is the Fermi distribution function. For 
a BCS superfluid gas with finite gap A the outcou- 
plcd density and current for each spin can be calculated 
using the BCS distribution function [53 /Bcs(£'n) = 

1 - en/v^JTA^tanh(yeITA^/2fcBT) where ^n = 

En — Ep, Ep is the Fermi energy and T the temperature 
of the system. 

The current is shown in Fig.[3for a spin-polarized nor- 
mal gas. The ID current in H69|) has an energy spread 
on the order of the Fermi energy Ep = Nfiuj^ due to 
the energy shifts of the contributions of different modes. 
The convolution with the transverse direction results in 
an even wider spread of the total current, which is caused 
by the tight confinement in the transverse direction, as 
explained in Sec. El It was predicted in |2J| and demon- 
strated experimentally |26l Isij that superfluidity in a 
Fermi gas leads to an energy shift in the outcouplcd cur- 
rent when a RF (or laser) field transfers atoms from one 
of the paired hyperfine states to another hyperfine state. 
The shift originates in the additional energy needed to 
outcouple the atoms that are forming Cooper pairs. In 
the BCS regime considered here the gap cannot exceed a 
small percentage of the Fermi energy. Such values would 
not create an appreciable shift in the current in Fig [3 
However, in a strongly interacting fermionic superfluid 
[sol I51II53.I53II the gap can reach values of the order of 
0.2Ep [261 Im | that would lead to an experimentally de- 
tectable shift in the current in Fig. [7\ 

The one particle correlation of the BEC was measured 
experimentally by outcoupling particles from the BEC 
with two different RF jloj. An equivalent process in a 
Fermi gas would lead to a density profile of the form 

n(r,t;z/i,i/2) - ^ /(S„)|V'out(ri, r; /ij/i -f S„)e'''"'* 

n 

+Vout(n,r;/ii/2 + S„)e*''=*p. (73) 
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been analyzed and a simple way to convolute the prop- 
agation along the gravitational and transverse direction 
has been determined. We have extended the formalism to 
calculate the current and outcoupled beam from excited 
modes. We have applied it to calculate the current and 
beam profile for a quasi- ID Fermi gas. We have shown 
that the interference pattern of atoms outcoupled with 
two different RF show oscillations that are not related 
to superfluidity as in the case of a Bose gas. The effect 
of fermionic superfluidity in the current and outcoupled 
beam has been discussed. 



APPENDIX A: RECURSION RELATIONS FOR 

THE CURRENT AND BEAM DENSITY FROM 

THE EXCITED MODES OF A TRAP 



FIG. 8: Density profile of the fermionic beam outcoupled 
with two different RF l-Bgrav — -Btrapl/ft and IkHz-f (_Egrav + 
Etr!!,p)/h. We show the density profile (eq. (I73II ) at the center 
of the beam for a Fermi gas with A'' — 70, u^ — 27rl.2 kHz, 
uj± — 27rl20 kHz, and pairing gap A — 3.36 kHz, 'y/h — 
lOOHz, rriK = 40 u. 



As demonstrated in Fig. |S1 the density profile at the cen- 
ter of the beam shows oscillations. For a Fermi gas, the 
oscillations are not related to superfluidity as in a BEC 
(see (123) and [l3|). In ID, the density profile lO is a 
sum of oscillatory functions in E+Fz (see eq. (j^ ). Each 
of the terms in the sum is shifted by huj,-, but because 
Ef «C Fz the period of oscillation is effectively the same. 
The 3D density profile H7U|) is a convolution of the ID 
density function with the energy of the transverse direc- 
tions that would result in oscillations on the scale E+Fz. 
Each term of the sum will contribute again by shifting 
the rescaled energy and as long as Ep <C Fz all of them 
will have effectively the same period of oscillation. As 
expected, the one particle correlation function does not 
show any effect of the superfluidity. The gap energy shift 
transforms into a shift 6z = A/F (see Fig. ISJ that only 
leads to a time shift in the density profiles. Fermionic 
superfluidity relies on the formation of atomic pairs and 
therefore one expects to see some effect of superfluidity 
only in the two particle correlation function 



VII. SUMMARY AND OUTLOOK 

A quantum mechanical theory of the atom laser based 
on propagator techniques has been presented and con- 
trasted with various approximation schemes. The anal- 
ysis of the current and outcoupled beam leads to a dis- 
tinction between small and big condensates. The experi- 
mentally observed structure and substructure [ij, [13 of 
the transverse beam proflle has been obtained using the 
T-matrix formalism, which includes the effect of interac- 
tions. 

The effect of a non-isotropic geometry of the trap has 



We derive eq. (|67|l for the outcoupled wavefunction 
from an excited state in a one dimensional harmonic 
trap. The generating function of the Hcrmite polyno- 
mials reads 



H„ - 

■.a 



gn 

—— exp(-t^ + 2tz/a) 



t=o 



Inserting it into Ea. (|66|l leads to 



i;l^,{n,z,E)^j I dz'G'^,{z,z\E)M^') 



= 7 



(Al) 



(A2) 



^ __ , dz'G(z,z',S)H„(-)e-- 



= 7 



V2"n! 



a" _t2jdz'G{z,z\E)e- 



r+2t^ 



Introducing a new variable z" = z' -\- 2ta and using the 
translation law for the Green function [lj|, Eq- (33), 
yields 



^^i^tin,z,E) 



(A3) 



7 



V'2"n! 

7 



V2"n! 



9" _t2^2^2jdz"G{z,z" + 2ta,E)e-^ 

m^"^ VS^i/4 

9" ^2jdz"G{z-2ta,z",E + 2taF)e~2. 



t=o 



'air 



1/4 



t=o 



The z" integration in the last expression yields the known 
outgoing wave function for n = 0, 

^l^t{n,z,E) 



7 



V2"n! 

7 



Qn 

dt 



-e' V'out(0,;2-2ta,^ + 2ta) 



V2"n! 

Qn 



4V2^5/4y;^;^g-8aV3+2a^eci(g_ 2|) X 



r-8a^t 



at" 



Ai(e - 4at) 



t=o 



B{n, z, E) [K{n, 2, -8a^, -2I3E + Aa^, -4a)] .(A4) 
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Here, wc separated the factors in front of the square 
brackets from the derivatives and used the K{. . .) no- 
tation of Ref. 1^, eq. (13): 



K{n,aL,a'L,'yL,5L) 



|Le^".*^+"l*Ai(7i + ,5it) 



t=o 
(A5) 
The expressions K(n) are readily calculated using a re- 
cursion relation derived in [53 that we show here for 
completeness. One can define 



K'{n,aL,a'i^,lL,5L) 



^e^"-*'+<*Ai'(7i4-M 



and simultaneously calculate K{'n) and K'{n): 



(A6) 



K{0) = Ai(7i), K{1) = a'iA-(O) + (5^/^1(0), 
K'{Q) - Ai'(7i), K'{\) = a'^K'iQ) + 5l1lK{0), 
K{n) = a'LK{n-l) + aL{n-l)K{n-2) 

K'{n) = a'^K'{n-l) + aL{n-l)K'{n-2) 

+5LlLK{n - 1) + 5l{n - l)K{n - 2). (A7) 



As pointed out in [55j the recursion method is unstable 
for \5l\ < 1/2 and for large 7^. This means that we 
cannot use this method to calculate the current and out- 
coupled beam profile for small trapping frequencies or 
large number of atoms. For small trapping frequencies, 
one could use the refiection approximation eq. H32(l . 
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